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Abstract 

An updated version of the integrated 
design program C/CARES (composite ceram- 
ics analysis and reliability evaluation of 
structures) was developed for the reliability 
evaluation of CMC laminated shell compo- 
nents. The algorithm is now split into two 
modules: a finite-element data interface pro- 
gram and a reliability evaluation algorithm. 
More flexibility is achieved, allowing for easy 
implementation with various finite-element 
programs. The interface program creates a 
neutral data base which is then read by the 
reliability module. This neutral data base 
concept allows easy data transfer between 
different computer systems. The new inter- 
face program from the finite-element code 
MARC also includes the option of using 
hybrid laminates (a combination of plies of 
different materials or different layups) and 
allows for variations in temperature fields 
throughout the component. In the current 
version of C/CARES, a subelement technique 
was implemented, enabling stress gradients 
within an element to be taken into account. 
The noninteractive reliability function is now 
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evaluated at each Gaussian integration point 
instead of using averaging techniques. As a 
result of the increased number of stress evalu- 
ation points, considerable improvements in the 
accuracy of reliability analyses have been 
realized. 

Introduction 

The development of propulsion systems 
for space vehicles as well as for hypersonic 
airplanes requires the use of engineered mate- 
rials such as laminated ceramic matrix com- 
posites (CMC’s). CMC material systems have 
low density and the advantage of maintaining 
strength at high temperatures, thus the 
demand for them in the design of high-temper- 
ature sections of aerospace structures. This 
demand requires focused research with specific 
attention paid to their inherent brittleness 
and variability in strength. Thus, analyses of 
components fabricated from ceramic materials 
require a departure from the usual determinis- 
tic design philosophy (i.e., the safety factor 
approach) prevalent in analyzing metallic 
structural components. Adopting the probabi- 
listic approach outlined herein allows the 
design engineer to account for variability in 
strength and decreasing bulk strength with 
increasing component size. 



As in any structural analysis of a com- 
ponent, the stress field must be characterized. 
Commercial finite-element programs (e.g., 
MSC/NASTRAN, MARC, ABAQUS, 

ANSYS, etc.) exist and have laminate analysis 
capabilities which allow the design engineer to 
determine the structural response of composite 
components subjected to thermomechanical 
loads. Coupling these finite-element programs 
to an integrated probabilistic composite design 
program that evaluates the structural reliabil- 
ity is an attractive analytical tool in compo- 
nent design. Based on the preliminary version 
of the test-bed computer program C/CARES 
(composite ceramics analysis and reliability 
evaluation of structures ), 1 which was coupled 
to MSC/NASTRAN laminate analyses, a 
next-generation algorithm was developed. 
Thus, access to additional finite-element pro- 
grams is allowed through the use of a format- 
ted neutral data base. In the current version, 
a subelement technique has been implemented. 
The reliability analysis is now performed at 
each Gaussian integration point instead of 
using averaging techniques. In a previous 
version of C/CARES , 1 element stresses were 
obtained by averaging the stress state at each 
node over the entire element. In finite- 
element analysis, nodal stress data are usually 
computed by extrapolating the stresses from 
the integration points; hence, the element 
stress tensor is approximate. This lack of 
accuracy becomes a difficulty in characterizing 
problems where high stress gradients exist 
over the element or where large element dis- 
tortion occurs. As a result of the new sub- 
element technique, the variation in stress over 
the element is taken into account, and con- 
siderable improvements in the accuracy of reli- 
ability analyses have been realized. 

In addition to a new subelement tech- 
nique, a modular program structure is intro- 
duced. Future changes and additions (e.g., 
implementation of more sophisticated failure 
models) are readily done without invoking 
large programming efforts. This version of 
C/CARES includes the option of using hybrid 
laminates, and it allows for variations in tem- 
perature fields throughout the component. 


In this report, the C/CARES interface for 
the finite-element package MARC is docu- 
mented. In MARC, several shell elements 
(which utilize different shell theories) can be 
used to model the linear as well as nonlinear 
structural response of a thin-walled compo- 
nent to thermomechanical loads. Most of the 
data required for the reliability analysis are 
supplied by MARC, including the ply stress 
states and temperatures at the integration 
points. The ply stresses are provided in the 
local material coordinate system. Improve- 
ments of the MARC finite-element data inter- 
face and the new version of the C/CARES 
program are evident when a comparison is 
made of the reliability results to previous 
MSC /N ASTRAN-C / CARES. 

Reliability Model 

A great deal of work has been done to 
develop theories of the micromechanics of fail- 
ure as they relate to ceramic matrix compos- 
ites . 2,3 Application of these theories directly 
to structures of reasonable size is not practical 
because of the large number of localized fail- 
ures and failure interactions which must be 
considered. For these reasons, the macro- 
scopic phenomenological theories that have 
been advocated 4 ' 9 rely on homogenization of 
material properties. However, the required 
strength parameters are now taken to be ran- 
dom variables. The reliability functions most 
often used were borrowed from existing deter- 
ministic models. These reliability functions 
can be classified as either interactive, in which 
case the stresses are assumed to act in a col- 
lective fashion when causing failure, or non- 
interactive, in which case all stresses are still 
considered, but individually. 

Noninteractive methods have the advan- 
tage of being evaluated in closed form and 
have been implemented in the past and cur- 
rent versions of C/CARES. Each unidirec- 
tional ply is considered independently and is 
allowed to possess five failure modes. (Delam- 
ination failure is not admitted in this 
approach.) In this model, the failure modes 
are characterized by tensile and compressive 
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strengths in the fiber direction (<r| ,max and 
^c,maxj; tensile an( j compressive strengths in 
the direction transverse to the fiber (< 72 ,max 
and <^ max ); and in-plane shear strength 
(a™?). These strengths are considered ran- 
dom variables that are characterized by 
three-parameter Weibull distribution func- 
tions. The reason for adopting this distribu- 
tion is the anticipation that laminated CMC 
materials will exhibit threshold behavior, at 
least in the fiber direction. Inclusion of the 
threshold stress parameter allows for zero 
probability of failure when the load is below a 
predetermined level. The determination of 
these strengths requires a large number of 
experimental tests. In order to establish a 
reliable data base for CMC material systems, 
further experimental efforts are necessary. 

To determine lower and upper reliability 
bounds, the weakest link theory (WLT, fail- 
ure of the weakest link causes failure of the 
whole component) and the strongest link 
theory (SLT, failure of the whole component 
will not occur until every link has failed) are 
applied, respectively. In both approaches, the 
different failure modes are assumed to act 
independently. If WLT is assumed, the reli- 
ability of the i th ply can be expressed by 

Rj = expj- V>i dvj (1) 


where V'i 13 the risk of rupture intensity (per 
unit volume), and V is the ply volume. 

The risk of rupture intensity ^ is a 
function of position x,y,z and of the local 
temperature T: 


^i(x,y^T) = 


<*i - 7 i(T)> 

^(T) 


K<T) 


+ 7i(T)> 




{ ^(T) j 


“ e i< T > 


( °2 ~ 7 2 (T)) 


4(t) 


«o( T ) 


(-*2 + 7 o (T)) 


£(T) 


« C ,(T)I 


f<kl*l - 7 i*(T)> 

*a(T) 


^«12(T)| 


(2) 


Note that 7 is the Weibull threshold param- 
eter, a is the Weibull shape parameter, and 
/3 is the Weibull scale parameter. In Eq. ( 2 ), 
these Weibull parameters are assumed to be 
temperature dependent. Superscript t repre- 
sents tensile, and superscript c denotes com- 
pressive conditions. Compressive stresses 
have a negative sign. In addition, o 1 is the 
ply stress in the fiber direction, o 2 13 the ply 
stress transverse to the fiber direction, and 
o 12 is the in-plane ply shear stress. The 
stresses are assumed to be constant over the 
ply thickness. Under the assumption of 
WLT, the component reliability is calculated 
by application of Eq. ( 1 ) to the whole 
component: 
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nply 

r wlt = n R i w 

i=l 

where nply is the number of plies in the 
laminate. Equation (3) represents a lower 
bound on the reliability. Thomas and 
Wet her hold 7 proposed an upper-bound reli- 
ability limit estimate by assuming SLT, 
where 

nply 

r slt = i-II I 1 - R i) (4) 

i=l 


Program Features 


Since the stress field throughout a com- 
ponent is required for reliability analysis using 
the C/CARES algorithm, the structural 
response of the component to thermomechan- 
ical loads must be determined first (via the 
application of finite-element methods). 

Because of their relatively small thickness in 
comparison to the overall dimensions of the 
structure, components fabricated from lami- 
nated CMC material systems are conveniently 
modeled by using shell elements. In the for- 
mulation of typical shell elements, classical 
lamination theory is usually adopted to 
describe their mechanical behavior. The 
finite-element results have to be filtered since 
only ply stresses and ply temperatures at the 
integration points are required for the subse- 
quent C/CARES reliability analysis. For con- 
venience, the program is split into two parts: 
a finite-element data interface program and a 
reliability evaluation algorithm. 


The finite-element data interface program 
interprets the data necessary for subsequent 
reliability analysis and translates the informa- 
tion into a formatted (ASCII) neutral data 
base. This concept enables easy data transfer 
between different computer systems. If the 
finite-element analysis and the reliability 
analysis are executed on the same computer 
system, the creation of an unformatted neutral 


data base is optional. In this case, faster and 
direct access to finite-element data is possible. 

The reliability evaluation algorithm per- 
forms computations necessary for determining 
structural reliability. In the subelement tech- 
nique, the risk of rupture intensity is evalua- 
ted at each Gaussian integration point instead 
of using averaging techniques. Stresses are 
exactly determined at the Gaussian integra- 
tion points where the local stiffness matrix has 
already been evaluated. The subelement tech- 
nique implies that each Gaussian integration 
point is referred to a subelement. The sub- 
element volume is defined as the contribution 
of the integration point to the element volume 
in the numerical integration procedure. Thus, 
each ply in the laminate is divided into sub- 
elements, and the risk of rupture intensity 
function rj) is evaluated for the corresponding 
stress tensor at the integration point. The 
number of subelements in each element 
depends on the order of integration chosen 
and on the element type. The increased num- 
ber of points for stress evaluation accounts for 
the variation in stress over the element. Since 
the stress state is assumed to be constant over 
the subelement, the volume associated with a 
stress evaluation point substantially influences 
the reliability evaluation terms (see Eq. (1)). 
Due to the subdiscretization in the subelement 
method, even problems with high stress gradi- 
ents or with distorted mesh geometries can be 
analyzed without remeshing the component 
(provided that the mesh density is sufficient 
to give accurate stresses). 

The reliability analysis requires the 
stress state of each ply subelement, the ply 
subvolume, the local temperature, and the 
temperature-dependent Weibull parameters. 
Based on this information, the subelement 
reliability is expressed as 

Rj Ub = exp |-V'j Ub Vj Ub | ( 5 ) 

where ^? ub is the value of the risk of rupture 
intensity within the subelement of the i th ply. 
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Similarly, V- ub is the volume of the sub- 
element of the i tb ply. 

The value of the risk of rupture intensity 
^? ub is determined by using Eq. ( 2 ). The 
Weibull parameters a, and 7 are calcu- 
lated with respect to the local temperature by 
linear interpolation of the material parameters 
specified in the control file. 

The reliability of each ply is determined 
as the product of all subelement reliabilities 
within the model: 


R i = n R f = ex p 

k 



dimension arrays) in each routine allows the 
user to easily adapt the program to the prob- 
lem size and to the available computational 
resources. Since the program has been devel- 
oped in a modular structure, future changes 
and additions (such as the introduction of 
more sophisticated reliability models) will not 
involve large programming efforts. 

The program requires two input files, a 
control file, and a neutral data base. Defined 
in the control file are the number of different 
materials, the temperature-dependent Weibull 
parameters, the maximum number of plies, 
and flags for printing and postprocessing (e.g., 
writing the element risk of ruptures to a 
PATRAN element file). The neutral data 
base is discussed in the next section. 


Thus, considerable improvements in the accu- 
racy of reliability analyses have been realized 
(the reader is referred to the numerical exam- 
ples that follow). After the reliability analysis 
is performed, the data for the risk of rupture 
intensity (a local measure of reliability) can be 
stored in a PATRAN element file. Thus, the 
design engineer is allowed to visualize critical 
regions in the structure by obtaining contour 
plots of the appropriate variables. The ele- 
ment risk of rupture intensity is approximated 
by averaging all subelements within the ele- 
ment and by using the following expression: 



where nsub is the number of subelements 
within one element {^f e is only used for 
PATRAN postprocessing). 

The current version of C/CARES is pro- 
grammed in standard FORTRAN 77 because 
it can be easily implemented on numerous 
computer systems from personal computers to 
mainframes. As a result of the subelement 
concept, a global structure of the algorithm 
has evolved. The introduction of one large, 
blank common array (instead of several fixed- 


Neutral Data Base Concept 

The formatted neutral data base allows 
C/CARES to be interfaced with several finite- 
element programs. At the present time the 
interface program to MARC has been imple- 
mented. Interfaces to MSC/NASTRAN and 
ABAQUS are being prepared. 

The structure of the neutral data base 
is optimized with respect to memory. The 
finite-element data are arranged within the 
neutral data base by using the following hier- 
archy: element groups, elements, and sub- 
elements. Element groups are already defined 
in the finite-element analysis. Elements with 
the same laminate structure are combined into 
one group. If reliability analysis is focused on 
a special part of a component (e.g., a region 
with stress concentrations), the corresponding 
elements can be combined into an element 
group. The element group data contain infor- 
mation regarding the number of elements 
within the group. Information pertaining to 
the elements includes the number of sub- 
elements (integration points) and the number 
of plies associated with each element. Finally, 
the subelement data contain information 
regarding ply identification, ply subvolume, 
ply stress state (assuming plane stress condi- 
tions), ply temperature, and ply material 
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type. Note that a subelement can have more 
than one ply. 

Typically, ply volumes are not included 
with standard finite-element output. Thus, 
the volume of each subelement (corresponding 
to a Gaussian integration point) is calculated 
by using the shape functions inherent to the 
element type. Currently, bilinear (four nodes) 
and biquadratic (eight nodes) element shape 
functions are utilized to determine the sub- 
element volumes. The subelement volumes 
are computed by using the following 
expression: 

^ply,sub “ ^ply,sub^ply ^ 

where A p ^ y gub is the midsurface area of the 
shell subelement (corresponding to the inte- 
gration point considered), and t ply is the ply 
thickness. 

The midsurface area of a shell element is 
calculated by integration in natural space (see 
Bathe 10 for details): 

A = J 1 J 1 detJ(r,s) dr ds (9) 

where J is the Jacobian operator, and r,s 
are natural coordinates. 

By applying a Gaussian integration 
scheme, this integral can be expressed as 

m m 

A = £ £ detJOfcspWjWj (1°) 

i=l i=l 

where m is the integration order, r.,Sj are 
the coordinates of the Gaussian integration 
points in natural space, and the W 5 s are 
weight functions. 

The Jacobian operator J is calculated by 
using the corresponding element shape func- 
tions, element connectivities, and global coor- 


dinates of the element nodes. The subelement 
area corresponding to the integration point 
rj,s- represents one summation term in 
Eq. ( 10 ). 

Numerical Examples 

To illustrate the improvements achieved 
by the updated version of C/CARES (coupled 
to the MARC finite-element program), results 
of a benchmark problem reliability analysis 
are compared to the previously completed 
analysis for the same problem where 
C/CARES was coupled to MSC/ 

NASTRAN . 11 Specifically, a laminated CMC 
plate with a central circular hole was sub- 
jected to prescribed edge displacements at 
room temperature (see Fig. 1 ). The elastic 
material properties are presented in Table 1 . 
The laminate consists of eight unidirectional 
orthotropic plies with a [0/d=45/90] g ply 
layup. Note that, because experimental data 
are lacking, assumed values for a, /?, and 7 
(which are temperature dependent) were used. 
In order to model this quasi-two-dimensional 
problem, appropriate degrees of freedom are 
retained in each analysis. 

Two mesh densities were used for the 
finite-element analysis of this problem. The 
first mesh contained 2520 elements; the second 
coarser mesh contained 1000 elements. A 
two-by-two integration scheme was used for 
both MARC analyses (four subelements for 
each shell element). The finite-element results 
were transformed into a neutral data base, 
and the reliability analysis was performed 
with the Weibull strength parameters given in 
Table 2 . In Table 3, the component reliabili- 
ties R WLT are given for both MARC-C/ 
CARES analyses and the original 
MSC/N ASTR AN-C /CARES analyses. The 
component reliability based on SLT is 
nearly unity, since the high reliabilities R| 
for plies 2 to 7 dominate the calculation of 
R S . In Table 4 the reliability Rj for each 
ply from the MARC analysis (dense mesh) is 
compared to the ply reliabilities from the 
MSC/N ASTRAN analysis. 
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The plies with a fiber orientation of 0° 
(plies 1 and 8) are subject to the highest 
stresses in the laminate. Consequently, these 
plies are the least reliable (see Table 4). In 
Figs. 2 and 3, the risk of rupture intensity 
e is depicted in the vicinity of the hole for 
ply 1 (fiber orientation 0°) and for ply 2 
(fiber orientation +45°), respectively. 

Because the risk of rupture intensity is 
directly dependent on the stress state, the 
stress concentration is clearly evident. 

The reliability analyses using MARC 
results yielded smaller component reliabilities 
for both the coarse and the fine mesh relative 
to the MSC/NASTRAN analyses. The 
decrease in reliability is a consequence of the 
more accurate reliability computations in the 
region of the stress gradient. Both MARC 
analyses yielded close values for the maximum 
stress in the fiber direction (within 1.6 per- 
cent); the value of c™*** 0 \ s 433.8 MPa for 
the coarse mesh, and the value of a ™ ax ’ d is 
440.8 MPa for the dense mesh. The maxi- 
mum stress in the fiber direction determined 
with MSC/NASTRAN, however, is smaller: 
the value of [ s 416.0 MPa (coarse 

mesh), and the value of is 423.0 MPa 

(dense mesh). These results may be attribut- 
ed to the averaging method used in 
MSC/NASTRAN to determine the centroidal 
stress for each element. Thus, the stress gra- 
dient in the component could only be approxi- 
mated. Stresses at the maximum stress 
location may be underestimated for the MSC/ 
NASTRAN analysis, resulting in a higher 
component reliability for this problem. The 
MARC results show that the mesh density is 
sufficiently fine, even for the coarse mesh, to 
determine stresses accurately at the integra- 
tion points. This ability to model a structural 
component with a coarse mesh by using sub- 
elements represents a considerable reduction 
of computational effort in component design. 

Conclusion 


Based on the previous version of 
C/CARES, an enhanced version of this inte- 


grated design program was developed for the 
reliability evaluation of ceramic matrix com- 
posite (CMC) laminated shell components. 

By dividing the original program into two 
separate modules, more flexibility is achieved, 
allowing for easy implementation with various 
finite-element programs. As a result of the 
introduction of a shell subelement technique, 
the number of reliability function evaluation 
points is increased. A numerical example has 
demonstrated improvements in the accuracy 
of reliability analyses. In order to achieve the 
same order of accuracy, the new version 
allows structural components to be modeled 
with coarser finite-element meshes than those 
required by the previous version of the pro- 
gram. The modular structure of the program 
allows the design engineer to incorporate addi- 
tional reliability functions without major 
programming efforts. Future versions of 
C/CARES will focus on this issue. 
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TABLE 1.— ELASTIC PROPERTIES 


Young’s modulus in fiber direction, E lf MPa 19 755.0 

Young’s modulus transverse to fiber direction, E 2 , MPa 1975.5 

Shear modulus, G lf G 1S , G 23 , MPa 700.0 

Poisson ratio, u 12 0.35 


TABLE 2.— WEIBULL STRENGTH PARAMETERS 


Failure mode 

Weibull strength parameter 


Shape, 

Scale, 

Threshold, 


a 

fi, 

7, 



MPa 

MPa 

Fiber direction 




Tensile, 

25 



Compressive, o-*» max 

35 



Transverse to fiber direction 




Tensile, crj ,max 

10 



Compressive, <Tj' max 

30 



In-plane shear, <x™ 2 X 

22 

420 

▼ 
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TABLE 3.— COMPARISON OF DIFFERENT 


RELIABILITY ANALYSES 


Analysis 

Number of 
elements 

Component 

reliability, 

R wkt 

MARC 

mm 

0.9005 

MSC/NASTRAN 


.9304 

MARC 


.8959 

MSC/NASTRAN 

— 

.9328 


TABLE 4.— PLY-BY-PLY RELIABILITY RESULTS 


' Ply 
number 

1 

1 

Ply reliability, 
2520 elements 

MSC/NASTRAN - 
previous C/CARES 

MARC - 
new C/ CARES 

1 


0.96587 

0.94652 

2 

+45 

.99999 

.99999 

3 

-45 



4 

mm 



5 

wsm 



6 

-45 



7 

+45 

▼ 

▼ 

8 

0 

.96587 

.94652 


/- 7- mm radius 

/ 



z l X 

Figure 1 .-Finite-element mesh for [0/±45/90] s laminate. Ply thickness t, 0.2 mm. 
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Risk of rupture 




Risk of rupture 
intensity, 

¥ 

.000000036 


.000000028 


.000000020 


.000000012 


.000000004 


Figure 3. — Risk of rupture intensity ^ around the hole for ply number 2 (+45° fiber orientation). 
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